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Abstract 

This work describes a novel image analysis approach to characterize the uniformity of objects 
in agglomerates by using the propagation of normal wavefronts. The problem of width uniformity 
is discussed and its importance for the characterization of composite structures normally found in 
physics and biology highlighted. The methodology involves identifying each cluster (i.e. connected 
component) of interest, which can correspond to objects or voids, and estimating the respective 
medial axes by using a recently proposed wavefront propagation approach, which is briefly reviewed. 
The distance values along such axes are identified and their mean and standard deviation values 
obtained. As illustrated with respect to synthetic and real objects (in vitro cultures of neuronal 
cells), the combined use of these two features provide a powerful description of the uniformity of the 
separation between the objects, presenting potential for several applications in material sciences 
and biology. 
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I. INTRODUCTION 



Structures defined by the agglomeration of basic elements are often characterized by the 
presence of objects and voids defined between the latter. Figure 1 illustrates such a situation 
with respect to neuronal cells grown in vitro on different substrata, namely fibronection (a) 
and tissue culture plastic (b). Several important properties of such structures can be inferred 
from the geometry of the objects and voids present in such images, which has motivated the 
growing use of imaging methods in physics and in biology (e.g. p, 0, 0, F° r example, 
the migration of nerve and glial cells and axon elongation in the developing nervous system 
lead to different patterns of cell agglomeration and separation ||. In material sciences, 
the geometrical properties of the voids have strong influence over the respective mechanical 
and electrical properties of the composite and can tell much about the composite formation 
process. The present work describes a new concept and methodology for characterizing the 
width uniformity of the objects or voids in such agglomerates. 





FIG. 1: Neuronal cells grown on different substrata, namely fibronectin (a) and tissue culture 
plastic. 



The first important point is to properly and accurately characterize what is meant by 
uniformity of objects or voids. A possible interpretation would regard how these elements 
distribute themselves along the composite space. We could be interested, for example, 
in quantifying their translational invariance, a problem that could be suitably tackled by 
using the lacunarity measure described in 0,0,0]. The alternative problem considered in 
the present work is to quantify the uniformity of the width of objects or voids. Figure 2 
illustrates such a situation with respect to two hypothetical structures chosen for didactical 
purposes. The agglomerate in (a) presents a single void (the black connected region) and 



2 



objects (in white) that have similar widths, while the case in (b) provides a non uniform 
collection of objects. The importance of quantifying the kind of width uniformity seen in 
(a) is particularly relevant because it may provide valuable information about the processes 
intrinsically related to the composite formation as well as about the interactions between the 
basic elements of the composite. Ideally, it would be interesting to have a single uniformity 
index that could be assigned to each agglomerate image. The current work proposes an 
effective solution to such a problem that can be extended also for analysis of connected 
components (i.e. objects or voids) involving connected holes such as the voids in the two 
situations in Figure 2. 




1(a) 



1(b) 



FIG. 2: Two types of agglomerates exhibiting different width uniformity of objects (represented in 
white). 



Firmly based on concepts from mathematics and physics, wavefront propagation ap- 
proaches have been extensively used for several purposes in image analysis and applications 



(e.g. 
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14j). Typically, such methods involve the numeric propagation, along the object 



contour normal, of a wavefront, starting from the object of interest, in such a way that shocks 
between portions or individual fronts define the medial axes of the objects under analysis. 
A shock is henceforth understood as the collision between two different points along the 
propagating wavefront (s). Medial axes are often defined as the set of points corresponding 
to the centers of circles maximally inscribed into the object of interest, which coincide with 
the shock points. Also known as skeletons, such axes can be informally said to correspond to 
the "middle" part of the objects (see Figure 5). Although alternative approaches to medial 
axes estimation (e.g. ^^) can be used, the current work considers the numeric approach 
reported in ( |2j, y, |4|), which also allows a spatial scale parameter controlling the degree 
of detail of the obtained axes. The basic principle underlying the method for quantifying 
width uniformity described in the current paper is the fact that the distance values along 



the medial axes provide a natural characterization of the width distribution of objects even 
when they are not straight or present ramifications and holes. The potential of the use of 
the distances along the medial axes is illustrated in Figure 3. Any object in this figure is 
characterized by the peculiar situation that the same distance value is found along the re- 
spective medial axes. The concept of width can therefore be extended from straight objects 
to a more general set of shapes, including those that are curved or present branches. In 
brief, all shapes in Figure 3(a) can be said to have the same width. This class of objects has 
been studied in under the name of polyballs. While such objects with constant width 
represent a somewhat limited class of possible objects, the concept of width can be imme- 
diately extended to more general shapes by considering the distance value along the object 
medial axes. Figure 3(b) illustrates such an approach with respect to a more generic object. 
The distances d along the respective medial axis represent by the dashed line, provide a 
natural means to characterize the widths of distribution the object. 




FIG. 3: Several objects characterized by the same distances along their respective medial axes(a). 
The generalization of the width of a generic object in terms of the distances along medial axes(b). 

Such distributions of distances along the medial axes are henceforth taken as an accu- 
rate indication about the object width. While the histogram of the distances along the 
medial axes of all objects (or voids) in the agglomerate could be taken as a more compre- 
hensive statistical model for characterizing the object width, we restrict our attention to its 
mean and standard deviation, represented henceforth as w and a w , as well as the correlation 
coefficient v between these two measurements. Given its adimensional nature, the latter 
provides a particularly interesting candidate for quantification of the uniformity of the ob- 
ject's width. As will be illustrated in this article, these three measures provide meaningful 
and comprehensive geometrical characterization of the objects or voids in agglomerates not 
only in terms of width uniformity, but also as far as the size of those objects is concerned. 



This article starts by presenting the methodology suggested to estimate the medial axes and 
distances, leading to the three width features w , a w , and v. Next, the potential of these 
measures to quantify width uniformity is illustrated with respect to synthetic and real data 
corresponding to agglomerates of neuronal cells grown over four different types of substrata. 



II. METHODOLOGY 



Given an image of an agglomerate characterized by objects over a background (or any 
two distinct phases) like that in Figure 2(b), the first step of the proposed methodology 
involves identifying all connected components for the object or voids. A connected com- 
ponent is understood as the set of all pixels that are connected to its 4 or 8 neighbors 
1 

4|. This work adopts 4- neighborhood. In case we are analyzing voids, a single connected 
component is obtained for the agglomerates in Figure 2. At the same time, a total of 10 
connected components are obtained in case we are interested in the objects. The extraction 
of the connected components can be done by using standard region growing or inundation 
algorithms (e.g. 0|). Let P be any of the points of the component to be extracted. Such 
methods search for the neighbors iVof Pthat have the same pixel value as P. The procedure 
is repeated for the neighbors ./V until no more neighbors with the same pixel value are found. 
Once the connected components have been obtained, we need to estimate the multiscale 
medial axes of each of them. This procedure is illustrated with respect to one of the void 
connected components from the agglomerate in Figure 2(b), shown in Figure 4(a). First, the 
border of the connected component is extracted in clockwise or counterclockwise fashion. 
Traditional algorithms for border extraction (e.g. ^|) are used in the present work. The 
output of such algorithms is a spatially quantified parametric curve defined by the list of 
ordered coordinates x(i),y(i) , where i=l,2,...,N is the curve parameter corresponding to the 
order the border elements are visited by the border extraction procedure. Figure 4(c) illus- 
trates a portion of the extracted border of the connected component in (a), parametrized 
by integer values corresponding to the sequence in which the border was followed. The arc 
length s(i) along such a curve is then accurately estimated by using the spectral method 
described in j4|. This method involves obtaining the Fourier series for the parametric 
curve defined by the extracted border and calculating the derivatives by using the property 
x <-> $5 _1 {Q ; (s)(j (J }, where S is the Fourier transform and G a is the Fourier transform of 
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the unit area Gaussian smoothing function G a = l/icryhr) exp{— 0.5(^) 2 } with standard 
deviation a. The arc length can be immediately obtained from the calculated derivatives 



as s(i) = ^2 y^ 2 (0 +y 2 {i) ■ Additional information about derivative estimation through 

fc=i n 
spectral methods can be found in |4| . The value of a is henceforth fixed at 7 pixels, which 

minimizes the effect of the spatially quantized nature of the curve - recall that in digital 
images such curves are represented over the orthogonal lattice - without excessively affecting 
(smoothing) the curve. The estimated arc length for the object in Figure 4(a) is given in 
(d). Having obtained the border of the connected component reparametrized in terms of 
its arc length, a wavefront is propagated from the border, unfolding with constant speed 
and normal orientation. The shocks implied by such a traveling wavefront as a consequence 
of curvature variations along the object geometry define the medial axis of the connected 
component under analysis. The methodology proposed in 0,0,0], and reviewed below, is 
used in this work for such a purpose. The basic idea of this method is to propagate cir- 
cular waves centered at each of the contour elements, according to the Huyghens principle 
of wave propagation. Constant radial speed is assumed. Given that the objects in digi- 
tal images are represented over the orthogonal lattice, which implies specific anisotropies, 
it is necessary to consider only those distances that are representable in such spaces. For 
instance, the first such distances are 0, 1, a/2,2, ... . Once such distances, which have been 
called exact distances, are pre-calculated (see 0) and stored jointly with the respective 
relative position vectors with respect to the origin, the algorithm consists in visiting each 
of the border elements for each subsequent exact distance, and updating the empty cells in- 
dexed by the respective relative positions with the current distance value or the respectively 
associated arc length. In the former case, the resulting image corresponds to the distance 
transform of the object, which associates to each surrounding pixel the smallest Euclidean 
distance to the object. In case the arc length is updated, we get a labeled image that can 
be further processed in order to obtain label differences (see 0, Q]) corresponding to the 
intensity of the shocks between different portions of the traveling wavefront. Figures 4(e) 
and (f) illustrate the distance transform and propagated labels obtained for the object in 
(a). 

By thresholding the difference image with different values T, it is possible to obtain a 
family of medial axes with varying degrees of detail. The higher the value of T, the more the 
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FIG. 4: Illustration of the multiscale medial axes estimation by using level sets. The original con- 
nected component(a), its extracted border(b), parameter values along a piece of such a border(c), 
arc-length along the parametrized contour(d), distance transform(e), and propagated labels(f). 

details are filtered out. The application of the above methodology for medial axes estimation 
to the problem of width uniformity involves selecting a suitable value of T so as to remove 
unimportant small scale detail such as vertices along the object. Figure 5 illustrates several 
medial axes obtained for varying values of T. The specific choice of the threshold value 
should take into account the resolution and specific demands imposed by each application. 



FIG. 5: Medial axes obtained for T = 5 (a), 10 (b), 20 (c) and 50 (d). 



If connected component presents holes, as is the case with the composites in Figure 1, 
each hole has to be identified by using border extraction algorithms and labeled individually 
with specific labels (i.e. a single distinct label is assigned to each border). By applying the 
wave propagation method to such borders, the generalized Voronoi diagram of the borders 
is obtained, and the multiscale medial axes are then calculated inside each of the Voronoi 
cells by U8ing the 8 a,ne m etho, desenbed above Q. The g ene r ahzed Voroaol d lagram 
corresponds to a partition of the image space in such a way that all points in each obtained 
region are closest to the respective object (one of the image connected component) than to 
any other object. Recall that the distance between a point and an object is defined as the 
smallest distance between the point and any of the points composing the object. Figure 
6 illustrates the generalized Voronoi diagram and the distance values along the borders 
(fixed-scale medial axes) determined by the Voronoi cells. 




FIG. 6: The generalized Voronoi diagram (b) for the composite in Figure 1(b) and the respective 
distance transform values (a) obtained along the Voronoi borders. 

At this point, it is possible to derive a distance weighted medial axis by assigning to 
each point of the selected medial axes the respective distance transform values. The set 



8 



of distance values found along the medial axis are henceforth called medial axis distances 
or skeleton distances. The methodology for width uniformity characterization involves the 
consideration of the distance values along all the medial axes in the agglomerate. The mean 
and standard deviation of the distances over every medial axis point of the object are then 
estimated. Table I gives these values, as well as the correlation coefficient, for the two images 
in Figure 1 considering objects and voids. 

TABLE I: The mean, standard deviation and correlation coefficient obtained for the voids and 
objects with respect to the two composites in Figure 1 





w 




V 


Fig. 1(a), objects 


2.06 


0.62 


0.30 


Fig. 1(a), voids 


6.57 


3.71 


0.57 


Fig. 1(b), objects 


4.37 


2.29 


0.52 


Fig. 1(b), voids 


6.09 


3.25 


0.54 



Several interesting facts result, or can be inferred from such measures. To start with, 
the mean width w obtained for the objects (cells) in Figure 2(a) is less than half of the 
value obtained for the objects in Figure 2(b), which is in full agreement with the visual 
interpretation of those images. At the same time, the mean widths obtained for the voids 
are similar in both images. The standard deviations obtained for the objects in those images 
are also in full agreement with the fact that the width dispersion of the objects in (a) is 
much smaller than that characterizing (b). Similar values are obtained for the voids. The 
coefficient of variation v provides an adimensional characterization of the observed width 
dispersion differences for the composites in Figure 2. 



III. APPLICATION TO NEURONAL AGGLOMERATES 

During the development of the nervous system the distribution and connections of cells is 
determined by a variety of factors, including the adhesion of cells and axons to extracellular 
matrix material such as laminin and fibronectin \^ little is know, however, about the 
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determinants of cell shape. Some aspects of cell-matrix adhesions and their relationship 
to cell shape can be studied in vitro by plating nerve cells derived from a single common 
precursor cell (clonal cells) on different surfaces and asking how the different surfaces regulate 
the morphology of the cell. Figure 7 shows that there are great differences in the shape of 
a clonal nerve cell line called B103 grown on polylysine, laminin, fibronection and tissue 
culture plastic. 400 X 319 pixel images of the cultured cells were obtained by using a Leitz 
DMIRB inverted microscope and an Openlab (Leica) image acquisition system connected 
to a personal computer. Cell shape can be assigned a numerical value and questions asked 
about how the shape of cells grown on the different substrata differ from each other. The 

n 

original gray- level images, such as those illustrated in Figure 1, were segmented |4J in order 
to separate the cells from the substrate, yielding binary images where the voids are marked 
as "0" and the cells (or objects) as "1". Attention is concentrated on the analysis of void 
uniformity. The segmentation was performed by automated segmentation (thresholding) 
followed by manual editing of the obtained structures in order to remove artifacts such as 
the presence of spurious particles and shadows. Figure 7 presents examples of binary images 
obtained by the segmentation methodology described above. Each of the connected voids 
were isolated by using a standard region growing approach |4j, and the multiscale medial 
axes were obtained by using the level set-based approach as described in Section 2. Specific 
medial axes were obtained by selecting the spatial scale to 10. 




20: w = 7.94 and a w = 3.62 32: w = 5.51 and a w = 2.67 
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56: w = 13.35 and a w = 5.26 



70: w= 10.09 and a w = 5.53 



FIG. 7: Binary images (after segmentation of the neuronal cells, shown in white, the voids are 
consequently black) showing examples of neuronal agglomerates grown over different substrata 
and respective values of w and cr w .(&): fibronectin; (b): laminin; (c): tissue culture plastic; (d): 
poly lysine. 

Figure 8 shows the phase space (or scatterplot) obtained by considering the mean and 
standard deviation of the void medial axes distances. The first important feature in this 
plot is the relatively strong overall correlation coefficient, namely 0.86, between the two 
considered measures. While an almost linear separation was obtained between classes 3 
(laminin, triangles) and 4 (polylysine, circles), the other two classes led to more dispersed, 
overlapped clusters. Regarding the mean value of w, calculated for each class, we have 
meaniw^) < meaniw?) < mean(w\) , indicating that the cells grown over these types of 
substrate tended to be progressively less packed. This is probably a reflection of a competi- 
tion between cell-cell and cell-substratum adhesion on each surface cells that adhere better 
to themselves than to the surface on which they are grown tend to aggregate or, in the case 
of neurons, form fascicles Q|. Therefore, cells grown on more adhesive surfaces such as 
tissue culture plastic would tend to be spread out on the surface, while those on less adhesive 
substrate would form clumps as observed in Figure 7. The correlation coefficient obtained 
for each of the four clusters are respectively 0.91, 0.73, 0.84 and 0.87, indicating substantially 
distinct uniformities characterizing each substrate. Whereas the voids in class 1 are char- 
acterized by the highest width uniformity, those in class 2 exhibit the most heterogeneous 
widths. 
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FIG. 8: Scatterplot showing the distribution of the neuronal agglomerates. The clusters obtained 
for class 3 and 4 are clearly linearly separable, while the other clusters are characterized by sub- 
stantial overlap. 

IV. CONCLUDING REMARKS 

This article has reported a novel and comprehensive approach to characterizing width 
uniformity in objects and voids typically found in images of agglomerates. The methodology 
is based on the objective characterization of "object width" in terms of the distances found 
along its medial axes, which are obtained by using a wave front propagation scheme. 
In this way, it becomes possible to speak objectively about the width of a wide class of 
shapes, including those that are curved or present branches. Such an approach to width 
characterization is further enhanced by the fact that most objects can be approximated 
or decomposed into varying radius polyballs. An accurate numerical approach has been 
adopted for the estimation of the medial axis distances that involves object isolation, border 
extraction, arc length reparametrization through a spectral approach, and the calculation 
of multiscale medial axes by using wavefront propagation. Several of these methods were 
only recently described in the literature and are applied in combined fashion for width 
estimation for the first time in the present work. The potential of the proposed framework 
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was illustrated and corroborated with respect to synthetic agglomerate images and real 
data concerning neuronal cells grown on different substrata. By examining a clonal nerve 
cell line grown on different protein substrata it was possible, using this new methodology, 
to quantify the morphological differences caused by the different substrata. It is shown that 
when nerve cells are grown on different surfaces, there is a dramatic correlation between the 
shape of the obtained voids and the tendency of the cells to associate between themselves. 
The use of the above methodology combined with complementary agglomerate measures 
such as the fractal dimension and lacunarity provide potential for a rich and comprehensive 
characterization of the geometry of the shapes and agglomerates of both biological and 
physical objects. 
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